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Abstract 

We propose a (time) multiscale method for the coarse-grained analysis of self- 
propelled particle models of swarms comprising a mixture of 'naive' and 'informed' in- 
dividuals, used to address questions related to collective motion and collective decision- 
making in animal groups. The method is based on projecting the particle configuration 
onto a single 'meta-particle' that consists of the group elongation and the mean group 
velocity and position. The collective states of the configuration can be associated with 
the transient and asymptotic transport properties of the random walk followed by the 
meta-particle. These properties can be accurately predicted by an advection-diffusion 
equation with memory (ADEM) whose parameters are obtained from a mean group ve- 
locity time series obtained from a single simulation run of the individual-based model. 

keywords continuous time random walks, anomalous transport, collective animal behav- 
ior, non-Markovian stochastic processes, self-propelled particle models. 

1 Introduction 

Self-propelled particle models (SPP's) are a class of agent-based simulations that have been 
used over the last three decades to explore questions related to various kinds of collective 
motion in animals, including insect swarms, bird flocks and fish schools [U EHl 1221 EEl HZl 
\W\ HSl mi- In these models, each individual in the (finite) population is represented 
by a particle that moves with constant speed in two or three-dimensional Euclidean space 
or a 2-dimensional torus. All particles update their orientations according to a set of lo- 
cal averages of the current state of the configuration. These local averages are simplified 
representations of individual behaviors that depend on 'social interactions' -avoidance of 
collisions, attraction, and orientation alignment- which result in the remarkable property 
of cohesive collective motion; i.e. the particles move about in space, yet they appear to 
move as a single object, resembling the motion of real flocks [H |50l [16]. Errors made by 
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the individuals as they estimate these quantities are modeled by a random rotation of the 
output of this averaging procedure. 

More recently, SPP models of flocking have been introduced in the context of collective 
decision-making to illuminate the question of how groups of agents achieve consensual deci- 
sions without the need of a central control [161 1131 113 [I2l SI] • Each of these decisions can 
be associated with a variety of collective states, which typically involve switching between 
mobile/immobile regimes [3S], rotation or milling J^E\, motion with a directional bias [12], 
or a combination of these [12] . A directional bias is relevant when critically important infor- 
mation, for instance the location of a resource, a predator or a migratory route, is available 
only to a fraction of the population [151 112] • 

[T^ explored this situation using a modified version of earlier models of swarming [H EOl 1221 
[TB] , where the main innovation consisted of dividing the population into two types. The first 
of these, called 'naive', follow only the social rules mentioned earlier (avoidance, attraction 
and alignment). The second kind, dubbed 'informed', also obey the social interactions of the 
naive individuals, but weigh the social output with an orientation bias along a single 'pre- 
ferred' direction, which in this study is identical for all informed individuals. This orientation 
bias can be regarded as a simple representation of access to privileged information. Collec- 
tive decision-making is understood in this context in terms of the ability of the informed 
sub-population to transfer their orientation bias to the whole group while simultaneously 
preserving group cohesion. 

Despite the recent explosion of SPP models in the literature, our understanding of these 
systems still remains limited. Central challenges are related to our ability to characterize 
efficiently and meaningfully the dynamics of each collective state, and critically, their depen- 
dence on the parameters of the individual-level model. We identify three distinct approaches 
to address this problem; namely Monte-Carlo simulation, continuum models, and 'hybrid' 
multi-scale approaches. 

The first (the Lagrangian approach) is mainly computational and consists of moving with 
each individual particle. Macroscopic summary statistics describing the various collective 
states are obtained from averages based on a large number of independent simulation runs, 
or a single time series when ergodicity is a reasonable assumption. These average quanti- 
ties usually include the mean group velocity [XT], the mean angular momentum [TU], mean 
switching times between mobile/immobile states [3^ or the 'accuracy' of the decision-making 
process [IHl [3H]. Other state variables of interest link collective states to geometrical prop- 
erties of the fiock, like the group elongation [15] or its aspect ratio [8] . 
The second method (the Eulerian approach) focuses on continuum models for the density 
and velocity fields. It has the advantage that in some cases analytical results linking the 
microscopic to the macroscopic can be rigorously derived. In addition to this, the numerical 
solution of the model for large or small population densities has the same computational cost, 
and the mechanisms that generate the collective patterns can often be clearly distinguished 
in the various terms in the model, which provides some degree of parsimony that approaches 
based solely on Monte-Carlo simulations cannot emulate. Continuum approximations have 
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been used to approximate discrete SPP models mainly to study collective motion that is not 
cohesive [17]; i.e. the population lives in a spatial arena with periodic or reflecting bound- 
aries but does not form a single distinct group. Instead, particles move about freely forming 
and dissolving groups of various sizes (i.e. fission-fusion dynamics), and collective motion 
is detected as a non-vanishing population average of the velocity. These continuum models 
are obtained through heuristic reasoning based on careful observation of system symmetries, 
or the invocation of conservation laws [IHl ISl ESI [IS]- 

Although substantial progress has been made with Eulerian (continuum) approaches, par- 
ticularly for swarming microbial populations [21 |5T], there are still a number of issues that 
preclude their widespread use. First, the use of heuristics does not clarify the dependence 
of the macroscopic parameters on the individual-level model. Although some continuum 
models have recently been derived formally from the individual-based model via a limiting 
process (usually large population size), the theoretical progress is made at the expense of 
great simplifications which restrict strongly their biological relevance. For instance [1] and 
[9] each derived continuum models in the limit of large population sizes, but restricted the 
individual-level interactions to a single type of social interaction, specified via a potential 
function [9j, or a velocity average [1]. Second, they usually require very large population sizes 
in order to be meaningful, which is problematic for models of flocking in groups involving 
tens or perhaps hundreds of individuals. In this situation the finiteness of the population 
size plays a fundamental role in observed transport properties (e.g. the group tends to move 
more slowly as the population size increases) [TH [71 [H [25], [531 [S3 • 

The third is the hybrid multiscale approach, which attempts to bridge Monte-Carlo simu- 
lations and continuum models. It is based on assuming the existence of a continuum model 
for some relevant coarse-grained state variable or 'reduction coordinate'; for instance Non- 
linear Advection-Diffusion Equations (NADE) with density-dependent coefficients [25], or 
Fokker-Planck equations with a non-linear potential [331 1201 HSl |35l [TTl [60] , which serves as 
a model template. The unknown fluxes and coefficients in the macroscopic template are es- 
timated from a computational experiment, which usually consists of a single -and relatively 
short- simulation run of the microscopic model. These estimated quantities are substituted 
into the unknown terms in the macroscopic model, which is then analyzed by means of the 
appropriate suite of classical continuum methods, numerical or analytical. 

In this study, we use this latter approach to explore the ability of Continuous Time Ran- 
dom Walks (CTRW) [lH [301 [lH [3] , its associated continuum counterpart, the Advection- 
Diffusion Equation with Memory (ADEM) -also known as the Generalized Master Equation 
(GME)- as a model template for the coarse-grained dynamics of cohesive collective motion 
and collective decision-making in self-propelled particle models of swarms comprising a mix- 
ture of individuals that have preferential access to critical information -the 'informed' type- 
and those who do not ( 'naive'). The ADEM generalizes the classical advection-diffusion 
equation to a non-local-in-time transport model via the introduction of a 'memory', a time 
weighting function proportional to the particle's two-time velocity autocorrelation function. 
The ADEM is a useful model of anomalous transport that arises when the underlying ran- 
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dom walk possesses a wide distribution of transition rates [301 EH [HI EHl E]- The multiscale 
method we propose is based on coarse-graining the full SPP configuration into a single 
'meta-particle', that consists of the group elongation (as a measure that the group remains 
cohesive) and the mean group velocity and position. The various types of collective states 
displayed by the group can then be related to the transport properties of the meta-particle's 
random walk, under the assumption that the pdf of the transition density for the meta- 
particle's position follows an ADEM. 

We illustrate the method for the case of a 2-dimensional SPP model introduced earlier by 
[15j for a single informed direction, but the approach is quite general in the sense that it can 
be applied to any individual-based model of movement for which the biologically meaningful 
coarse variables are the mean group position and velocity, and that the effective distribution 
of jump lengths at each transition event has finite moments of all orders. The multiscale 
approach for collective motion based on the ADEM complements /oca/-in-time multiscale 
approaches for a similar class of individual-level models explored earlier [35| |25| [60] . For 
instance, the ADEM can predict correctly the transport properties even when the individual- 
level model has a strong alignment rule, which is precisely the main limitation of the otherwise 
successful method based on non-linear advection-diffusion equations {23] . This results from 
temporal correlations in velocity fluctuations induced by the alignment rule that persist over 
macroscopically relevant time scales, a property that can not be captured by local-in-time 
Markovian models, but can be dealt with via the introduction of a memory term. 
CTRW theory generalizes the classical Random Walk (RW) as a microscopic model better 
suited for problems in anomalous transport, which is usually detected when the mean squared 
displacement (msd) does not scale linearly with time over a wide range of time scales. The 
anomalous properties can frequently be attributed to the presence of a wide distribution of 
transition rates (or also in the jump lengths), which leads to persistent temporal correlations 
in velocity fluctuations. It is the presence of time correlations in velocity that ultimately 
leads to anomalous transport [HI [301 [3H [SH]- The variability in transition rates can be 
attributed in real systems to spatial disorder in the medium, as is the case in tracer trans- 
port in porous media |3]. The presence of spatial disorder in the medium creates localized 
structures that can trap the particle for long periods of time, or force it to move ballistically 
by confining its motion along a corridor. The resulting particle motion consists of alternat- 
ing bursts of ballistic motion, apparent brownian motion, and a stagnant phase where the 
particle moves very slowly, if at all. This resembles the dynamics of the group meta-particle 
in SPP models of swarms, which typically consists of bursts of alignment in the particle 
orientations that lead to advective flights at the group level (the slip phase), alternating 
with regimes of slow motion when the particles lose their alignment and the mean group 
velocity drops sharply (the stick phase). The power of CTRW [HI [3] and effective medium 
theories of random motion in disordered media [31], lies in that the spatial inhomogeneities 
in the medium responsible for the anomalies in transport properties are not modeled explic- 
itly. Instead, their effect is summarized statistically in terms of the effective distributions of 
jump lengths and waiting times that define the random walk. The key innovation of CTRW 
theory is that the random walk does not proceed by fixed spatial and temporal increments. 
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but these become instead random variables, defined by two probability densities, which are 
usually assumed independent in applications. The first is the distribution of jumps in space 
A(^), which prescribes the length of the jumps between locations at each transition event. 
The second is a clock that regulates the times elapsed between transitions, known as the 
distribution of waiting times iP{t). A thorough discussion of modern CTRW theory and its 
role in models of anomalous transport can be found in a recent review by [32] • 
It can be shown [nil SH EDI ESI E] that when the distribution of jump lengths can be expanded 
in a Taylor series and the distribution of waiting times is an arbitrary probability density 
function, the transition probability density p(x, t|0,0) for finding a particle around position 
X G at time t given that it started at the origin at time zero, obeys a modified version of 
the advection-diffusion equation that is non-local in time, known as the Advection-Diffusion 
Equation with Memory (ADEM) 



- [ M{t- s) [va ■ Vp(x, s) - Da : V Vj9(x, s)] ds 
Jo 



dpjx, t) 
dt 

p(x,0+) = 5(x), xGM^ tG 



where va is the effective drift vector. Da the diffusivity tensor, and the colon operator is the 
inner tensor product 

A : B = Trace{B^ . A} = ^ A^jBij. 

hj 

The transport coefficients va and Da are determined respectively by the ratio of the first 
two moments of the jump distribution to the mean of the waiting time distribution, or the 
median when il){t) does not have a finite mean [3]. The memory function M{t) has two 
equivalent interpretations; it is closely related to the distribution of waiting times jH], but 
it can also be shown to be proportional to the velocity time auto-correlation function of the 
moving particle [28l |32l [58] . 

E[vi{0)vi{r)] = 2DiM{T), 

where vi is the velocity along the xi direction, Di is the diffusivity along xi and M(t) is the 
memory kernel that prescribes the decay of correlations (see Section |4]for additional details). 
This model constitutes the basis for effective medium theories of anomalous transport in dis- 
ordered media, where the spatial disorder in the medium is replaced by an ordered model 
with memory of the form Q [2HI EZl E21 El EI] • 

The stochastic dynamics of the meta-particle associated with the SPP model of flocking ex- 
plored here has a striking resemblance to that which motivated the development of the theory 
for anomalous transport in heterogenous media based on the CTRW and the ADEM. In SPP 
models, the wide range of variability in transition rates cannot be attributed to spatial disor- 
der in the medium, but arises instead from stochasticity in the alternating (slip/stick) types 
of collective behavior. Even though the source of variability is quite different, this does not 
seem to matter provided that transport can be modeled in terms of effective distributions 
of jump lengths and waiting times and their associated ADEM. Our goal is to exploit this 
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analogy to propose an ADEM as a continuum 'model template' for the dynamics of the 
position pdf of a swarm centroid. The functional form of the memory and the transport 
parameters in the ADEM template can be estimated from a single mean group velocity time 
series obtained from a simulation run of the SPP model. The resulting fitted model can 
be then used to explore the dependence of the collective behaviors on the parameters that 
determine the individual-level model, particularly the strength of the bias of the informed 
sub-population, the total population size, and the proportion of informed individuals. 

In the ADEM approach, the memory is the fundamental object that encodes all the trans- 
port coefficients, the various transport regimes and their characteristic timescales [361 ES] . 
When the spatial distribution of the disorder is known |ST] or the Hamiltonian of the micro- 
scopic model [22] , it is possible to derive the memory in ([T]) from the microscopic dynamics. 
In general, one has to resort to simulations or experiments and subsequent function fittings, 
in order to obtain the velocity auto-correlation function. The non-linearities involved in the 
definition of the SPP seem to preclude the derivation of the velocity time auto-correlation 
function rigorously from the microscopic swarm model. We find from simulations that the 
memory kernel along the informed direction for SPP models can be very well fitted by two 
closely related functions. The first corresponds to a Gamma density, 

M{t) = ,^^—^t-^eM-t/ra), (2) 

which works well in swarms where there are no informed individuals present, but also when 
the proportion of informed individuals is small (and relatively low values of the coupling 
strength). The initial power law decay in (|2| leads to a sub-ballistic, super-diffusive tran- 
sient detectable in the mean-squared displacement. This power law behavior has an expo- 
nential truncation at a characteristic time scale that establishes the onset of the asymp- 
totic regime, which is dominated by diffusion in swarms with no informed individuals and 
a mixture of diffusion and advection (with constant drift) for groups that include informed 
individuals. We also find that the diffusion coefficient decreases with group size, and the 
time scale (tq) that determines the onset of the asymptotic regime increases with group size. 
The second, 'richer' situation, arises in informed swarms for high values of the bias along the 
informed direction, where the early time super-diffusive transient is followed by a regime 
where correlations oscillate before reaching the asymptotic state, which is also classical 
advection-diffusion. This additional regime requires a modification of the memory kernel (g 
in order to capture these oscillations. We find that a Mittag-Leffler function Ea,(5{z) with 
an exponential truncation [^ 159] , 

M{t) = t^-' E^^p [- {t/T.T] exp(-t/r,), (3) 

provides an excellent fit in this regime, at the cost of introducing two additional parameters 
(the exponent (3 and the time scale r^). We used these estimates together with the ADEM 
model in order to predict the behavior of the mean squared displacement (msd), i.e. the 
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second moment of the mean group position, which can be used to characterize the various 
types of collective behaviors and their characteristic time scales in terms of their effect on 
the meta-particle's transport properties. The functional forms themselves do not seem to 
change with group size, but only the parameters do. 

For the region of parameters where the group remains cohesive, we observed that there 
are two types of collective behavior that are shared by both naive (no informed individuals 
present) and informed groups. First, there is an anomalous super-diffusive transient at early 
times (the scaling exponent in the mean squared displacement lies between one and two) 
due to the prevalence of slip/stick dynamics over that domain of time scales. Asymptoti- 
cally, the msd scales linearly with time for naive groups (diffusion-dominated), but shows 
a sharp transition to quadratic scaling (advection-dominated) for informed ones along the 
informed direction, which indicates that on average, informed swarms diffuse, but also move 
with constant velocity over the longer time scales. This transition from linear to quadratic 
scaling allows the detection of the time scale at which the informed sub-population manages 
to transfer its orientation bias to the whole group; this time scale, or time to consensus, is a 
natural measure of the efficiency of the decision-making process. The magnitude of the drift, 
which depends on the degree of polarization of the particle orientations along the informed 
direction, is a straightforward macroscopic parameter for the degree of consensus. We also 
note that as the group size gets larger, the drift gets smaller for the same proportion of 
informed individuals and informed bias strength. Finally, the diffusion coefficient along the 
informed direction can be interpreted as a measure of the precision of the collective decision- 
making process -since it is a measure of the spread of an ensemble of swarm meta-particles- 
when compared with that of naive configurations. 

The resulting ADEM fitted from swarm simulation time-series is self-consistent in the sense 
that transport parameters estimated from the memory via a Kubo-Green relationship [36|l2l] 
coincide with those estimated from the moments of the jump and waiting time pdf 's of the 
associated CTRW for the three group sizes explored (A^ = 10,50,100), proportions of in- 
formed individuals, and strength of the bias along the preferred direction. We also discuss 
the phase diagrams for the transport coefficients estimated from this method, where we no- 
tice velocity-precision trade-offs: as the total group size gets larger, the decision-making 
becomes more precise at the expense of a slower mean group velocity. We also note that the 
time scale to consensus is invariant with respect to group size, and depends only on the pro- 
portion of informed individuals and the strength of the coupling along the informed direction. 

The paper is organized as follows: Section [2] introduces a slightly modified version of 
the SPP model with informed individuals of [15j, where we removed the constraint on the 
maximum turning angle that an individual can make during a time step. We then define 
the set of coarse-grained variables of interest, namely the group elongation, the mean group 
position, and the mean group velocity which we called the meta-particle. Simulation results 
are also shown, focusing on the mean squared displacement of the meta-particle as well as 
kernel density estimates of the probabilities of mean group speeds and orientations, finalizing 
with group elongation time series that detect when the group splits appart. These results 
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are later used to define macroscopic measures of collective motion and collective decision- 
making in terms of the transport regimes shown in the msd. Section |3]briefiy reviews known 
results from the theory of continuous time random walks (CTRW) ISBEQ], and its relation- 
ship to the advection-diffusion equation with memoy (ADEM) |3] that we use later as the 
macroscopic transport model for the transition density of the mean group position. Section 
|4] assumes that the random walk followed by the group meta-particle evolves according to a 
CTRW, and discusses the procedure used to estimate the memory and the transport coef- 
ficients of the associated ADEM, from a single velocity time series obtained from a run of 
the individual-based model. We compare mean squared displacements obtained from ensem- 
ble averages over simulation runs with those predicted by the fitted ADEM for which show 
analytical results for the time to consensus. The method is used to carry out a systematic 
exploration of the dependence of the macroscopic parameters -the diffusivity, the drift and 
the time to consensus- on the microscopic ones of immediate biological relevance; namely the 
relative proportion of informed individuals, the coupling strength, and the total population 
size. Some final remarks are presented in Section [5) 



2 Self-propelled particle model (SPP) with informed 
individuals 

Consider a population of j = 1, . . . , particles with positions Xj (t) in 2-dimensional Eu- 
clidean space. Each particle j moves with constant speed s along its orientation angle 6j{t) 
in [— vr, tt). We summarize this information as the (complex) particle velocity 

zj{t) =se''^^'\ 

The state of the population at (discrete) time t is represented by the configuration $t(A) 

<|.,(A) = {[x,(t),^,(t)]}, (4) 

where A is the region of observation. At each tick of the clock, the positions and orientations 
of each particle are updated according to, 

x,(t + At) = ^.(t) + s^^^^)At (5) 
Ojit + At) = (<l>t(D,)) exp(«AQ) 

where At is the time increment and ($j(Dj)) is a local average of the configuration restricted 
to an interaction region Dj centered around the j-th particle. The details of the averaging 
procedure are described in the collective motion rule below ( Figure [T]). Errors made by the 
individuals in their estimates of the local state of the configuration are modeled by rotating 
the updated orientation obtained from the local average by a random angle AQ, drawn from 
the wrapped Gaussian on the unit circle A/'u)(0, o"^ At) with mean zero and variance cr^At. 
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Figure 1: Interaction zones for a focal individual (blue) Xj (dot) with velocity Zj (arrow). 
The dots and the arrows represent other particles in the configuration (black). The region of 
avoidance is the interior of the circle of radius r^^,. The particles contributing to the region 
of alignment and attraction lie within the annulus of external radius Tat and internal radius 



The local average in (|5]) comprises two groups of rules. The first is based on the classical 
social interactions for collective motion, with parameters restricted to the domain in which 
the full configuration moves cohesively as a single object [U EDI [13 ED [IS]- The second is a 
steering rule proposed by Couzin et al [15j, that attempts to lead the motion of the group 
along a preferred orientation (3. This additional rule is followed only by a sub-population 
of 'informed individuals'. Whereas individuals that are not informed (called 'naive') update 
their orientations exclusively from the output of the social rules, informed individuals up- 
date their orientations according to a weighted average of the social interactions with the 
preferred direction. The weight of the bias along the preferred direction relative to the social 
rules is given by a 'coupling constant' cu, which is interpreted as a simple parameterization 
of an 'internal state' of the informed individual (e.g. starvation, detection of a predator or 
a resource). Collective decision-making is then understood in terms of the ability of the 
informed sub-population to transfer their orientation bias to the whole group. 

The three social interactions are: 1) avoidance of collisions, 2) attraction (centering), and 
3) alignment (polarization). Whereas the collective motion interactions are followed by all 
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particles, the steering rule is followed only by the informed sub-population of A^^ < 
particles, whose indices Jp = {ji, ■ ■ ■ , Jat^} are chosen uniformly from the set of indices of all 
the particles in the configuration J$ = {1, . . . , A^}. Both the number of informed particles 
as well as their indices remain fixed for all times once chosen at time zero. The particles 
that are not in the informed sub-group are called called 'naive'. Following Couzin et al [15] 
we have 



1. Collective motion rule 



(a) Avoidance of collisions 

We define the neighborhood of avoidance of the j-th particle Avj = B{rav,^j(t)) 
as the circular domain of radius Vav centered at Xj(t) (see Figure [T|. If the 
configuration restricted to the window Avj is not empty, the avoidance rule takes 
precedence over the other interactions. The avoidance rule prevents collisions by 
pointing the focal particle in the opposite direction of the centroid of the locations 
of the particles found within Avj, relative to the location of the focal particle Xj(t). 
The number of neighbors of the j-th particle in ^t{Avj) is 

N 

where the focal individual j is excluded from the count and /b(x) stands for the 
indicator function of some 2-D domain B, 

J / N f 1, if X G 5 
1^ otherwise. 

The vector pointing in the direction opposite to the centroid of the particles in 
Ay. is 

1 ^ 

The updated orientation due to avoidance is 

ej{t + At) = arg(dj)+Ag, (8) 
where AQ is a random angle drawn from A/'«,(0, cr^ At). 



(b) Attraction and Alignment 

If the configuration restricted to Avj is empty, we proceed to evaluate the align- 
ment and attraction updating rules. The neighborhood of attraction/alignment 
of the j-th particle is Atj = B{rat,^j(t) ), where B{rat,^j(t) ) is the circular do- 
main of radius Vat centered at Xj(t). The social interaction in this case is given 
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by the normalized vector sum over the positions (which determines the local at- 
traction vector), and the velocities, (which dictates the local alignment vector) of 
the neighbors. The number of neighbors in Atj is 



N 



^Ai, =5^/At,(Xfc(t)) 



k=l 



The contribution due to attraction is given by the vector d| pointing in the 
direction of the centroid of the positions of the neighbors relative to the focal 
individual 



N 



d|(t + At) = — ^ Im, (xfc(t)) [xfc(t) - x,(t) ] (9) 



and the contribution due to the alignment behavior comes from the average 
orientation of all the particles in Atj 



N 



d'j{t + At) = J2lAt,{Mt))zk{t). (10) 



fc=l 



The total contribution of the social rules is given by the vector sum of the nor- 
malized vectors associated with the attraction (|9| and alignment contributions 



(10) 



S(t + At) d%t + At) 

Mt + At) = -f ^ + ^^7 Aji- (11) 

' ||d(t + At)|| dp + At 



These two contributions are equally weighted in (11) but could be generalized 
so as to have different weights. In what follows we explore the former, mainly 
to explore the potential of the ADEM to predict the macroscopic dynamics in 
the presence of strong alignment, which has been shown to be problematic for 
Markovian models [25]. The updated orientation is given by the argument of the 
social interactions dj after rotating it by a small random angle AQ drawn as well 
from the wrapped Gaussian 

Ojit + At) = arg(d,)+Ag. (12) 

2. Steering rule for the informed sub-population 

If the index of the focal particle is in the list of informed indices J^, the updated 



direction is given by a compromise between the output of the social rules (12) and the 
informed individual's preference to move along the informed direction /3. This is given 
by the weighted vector average of these two contributions ( |2]) 

d*At + At) = Mj{t + At)+uj b, (13) 
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Figure 2: Updating rule for an informed particle. The updated direction corresponds to the 
normalized vector sum dj of the preferential direction vector b = (cos(/3), sin(/3)) rescaled 
by a factor u, with the unit vector pointing in the direction of the output of the social rules 
Uj = {cos{9j),sm{9j)). 



where a; is a weighting constant, uj is the unit orientation vector arising from the social 
rules (12) and b is the unit vector associated with the preferred orientation /3 and dj. 
The updated orientation is 

0,-(t + At) = arg[d*]. (14) 

Once all the particle's orientations are computed according to these social rules, the 
positions are updated according to (|5]). A summary of the parameters in the SPP 
model, together with the values used for the simulations are shown in Table [T] 



2.1 Simulation results 

Our coarse-grained analysis of the individual-based model consists of projecting the full 
configuration ^ onto a set of summary statistics that we dubbed the 'meta-particle'. We 
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Table 1: SPP model parameters 



Parameter 


symbol 


value 


units 


Radius of avoidance 




1.0 


m. 


Region of avoidance 


Av 


— 


m^ 


Radius of attraction 




5.0 


m. 


Region of attraction/alignment 


At 


- 


m^ 


Particle speed 


s 


1.0 


m/sec. 


Perception error 


a 


0.1 


radians 


Time step 


At 


0.1 


seconds 


Total population size 


N 


10 - 100 


individuals 


Informed population size 


Np 


-100 


individuals 


Coupling constant 


u 


0.0 - 0.6 


dimensionless 


Informed orientation angle 


/3 


0.0 


radians 



associate the stochastic properties of the meta-particle random walk to the various collec- 
tive states of the full configuration. For collective decision-making, we found that a useful 
projection consists of three state variables, the group elongation A(t), the mean group ve- 
locity v(t), and the mean group position x(t). The introduction of the group elongation is 
necessary in order to detect situations where the informed individuals leave the main group, 
something that occurs at high values of the coupling constant. In this situation collective 
decision-making is not consensual, since the informed individuals fail to lead the complete 
group along the informed orientation. 



2.1.1 Projections of the configuration 

The elongation is defined as the maximum of the set of two-point distances among all the 
positions of the particles in the configuration, 

A(t) :=max| ||x,(t) -Xfc(t)|| I Vj,A;GJ$}. (15) 

We restrict our study to values of the coupling constant u that preserve cohesive collective 
motion, in the sense that the whole configuration moves as a single entity [501 [TB] . This is 
tantamount to requiring A(t) to have a constant upper bound C, 

A(t) < C « oo. (16) 



If the property (16) is preserved, measures for consensual collective motion and decision- 
making can be developed in terms of the stochastic properties of the two other projections of 
the configuration which together with the elongation, define the configuration 'meta-particle' 
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Figure 3: Naive vs. Informed elongation dynamics for = 10 and 1 x 10® time steps. In 



both panels the blue graph corresponds to an elongation time series ( 15 ) from a configuration 



with no informed individuals (called a naive group). In both panels, the red graph shows 
the result of introducing a single informed individual, where u = 0.6. In the right panel, 
the black graph also corresponds to a configuration with a single informed particle, but for 
a higher value of the coupling constant u = 0.73. In both panels the group elongation A(t) 
remains bounded for all observed times for the naive configuration and the mild coupling 
{u = 0.3), indicating a configuration that moves cohesively. However, further increasing the 
coupling strength (black graph, u = 0.73) causes the group to split, as evidenced by an 
elongation that grows without bound. 



<^, = {A(t),x(t), v(t)}, (17) 

where the second element in the triplet is the mean group position, or configuration centroid 
x(t) 

1 ^ 

k=l 

and the third one is the mean group velocity \{t) or group polarization 

1 ^ 



In our simulations we observed that there is a non-trivial region of parameter space that 
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preserves cohesive collective motion (16) for both naive (no informed individuals present, 
Np = 0) and informed (at least one informed individual present, > 1) configurations. 
Figure [3] shows two scenarios for the dynamics of the elongation A{t). The left panel shows 
the situation where the cohesive collective motion property is preserved for two swarms 
of the same total population size (A^ = 10). Blue shows the elongation associated with 
a naive configuration, and red shows a configuration that includes an informed particle 
{u = 0.6, Nj3 = 1), observed during 1 x 10^ time steps. We see that the elongation associated 
with the configuration involving an informed individual tends to take higher values than in 
the naive case, but remains bounded. The right panel shows the effect of further increasing 
the coupling constant, {u = 0.73, black line) where the elongation remains bounded for some 
time (about 1 x 10^ time steps) after which it starts to increase, signaling that the group 
has broken apart. 

Figure |4] shows typical sample paths (blue lines) for the configuration centroid x(t) for a 
group of ten individuals and T = 1 x 10^ time steps. The full configuration at the end of 
the simulation is shown in the insets at the center of each panel, where red dots represent 
the locations of naive individuals, and blue the informed ones. Panel (a) corresponds to a 
configuration involving only naive individuals (iVg = 0), Panel (b) has one informed individ- 
ual {Np = 1) with a coupling constant u = 0.6. Finally, the lower panel (c) shows the case 
where five informed individuals are present (A^/3 = 5,0; = 0.6). The insets to the right show 
the same sample path over smaller spatial scales. In the inset of panel (a) we se evidence of 
separate clusters -where the group moves very slowly due to a lack of polarization (the slip 
phase)- connected by advective flights due to bursts of phase alignment (the stick phase). 
This behavior signals the slip/stick dynamics characteristic of these systems [35], and resem- 
bles the behavior of tracer transport in porous media with a preferential flow direction p], 
where the corridors that confine the tracer play a roughly similar role to the polarization 
bursts in SPP models that lead to ballistic flights, alternating with traps that slow the tracer 
motion, akin to the slip phase in the SPP. Panel (b) shows the result of adding one informed 
individual (blue dot) with a relatively high value of the coupling constant {u = 0.6) where 
without loss of generality we identified the preferred direction with the positive xi axis. The 
introduction of a single individual is enough to break the orientation symmetry of the naive 
case, resulting in a motion bias along the preferred direction, and an disentanglement of the 
clusters that appear in the naive case. Adding more individuals for the same coupling con- 
stant leads to a higher mean velocity. In addition to this, the motion develops an oscillatory 
behavior along the coordinate perpendicular to the informed direction. 

Figure [5] shows kernel density estimates of the probability density of realized mean group 
velocities v obtained from a single time series (T = 3 x 10^ time steps) collected after a 
transient of 1 x 10^ time steps, for swarms of ten individuals and three different informed 
regimes. The left panel corresponds to the probability density of mean group orientations 
arg(v), and the right panel to the modulus (or mean group speeds) |v|. We see that in the 
naive swarm (blue graph in the left panel), mean group orientations are chosen uniformly 
from [— vr, tt) at all times. However, this rotational symmetry is broken upon the introduc- 
tion of a single informed individual, which yields a symmetric density centered around the 
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Figure 4: Group centroid sample paths. The blue rugged line in the three panels corresponds 
to an individual sample path of the configuration centroid x(t) up to t/ = 1 x 10^ time 
steps starting near the origin. The parameters shared in all three cases are the population 
size = 10, as well as the collective motion parameters, given hj a = 0.1 radians, r^j = 
5.0m., Tav = 1.0m., s = 1.0m. /sec. At = 0.1 sees. In (a) all the individuals are naive, in 
(b) there is one informed individual with coupling constant u = 0.6, and in (c) there are 5 
informed individuals, also with u = 0.6. The preferred direction is the positive xi axis. 
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Figure 5: Empirical argument (left) and modulus (right) pdf 's for the mean group velocity 
(19) for a configuration of ten individuals. Both panels show kernel density estimates from 
a single velocity time series of 3 x 10^ data points. The blue graph corresponds to the naive 
configuration (A^^ = 0), red to a group with one informed individual (A^^ = l,uj = 0.3), and 
black shows the results for a configuration involving five informed individuals, and the same 
value of the coupling constant (A^/? = 5,co = 0.3). 



informed direction (red, left panel). The existence of peaks reflects a tension between the 
slip/stick dynamics and the biased motion along the informed direction p. Overall, the 
group moves along the informed direction, but slip/stick bursts are strong enough to par- 
tially counter that bias by trying to recover the rotational symmetry. Further increasing the 
number of individuals (black, left panel) leads to a unimodal density concentrated around 
the informed direction. Kernel density estimates for the mean group speeds shown in the 
right panel v show comparatively less variability between the naive and informed regimes. 
This is to be expected, since the mean group speeds arises mainly from the social rules, and 
the steering rule is designed to introduce an orientation bias, but has little effect on the 
modulus. There is however a tendency to move with higher speeds as informed individuals 
are introduced. Critically, the range of variability in speeds is very wide and practically 
covers the full range of possible values (the individual particle speed is s = 1 m/s, which 
constitutes an upper bound for the mean group speed). 



2.1.2 Coarse variables and collective behaviors 

Measures for the collective behaviors that are macroscopically relevant can be defined in 
terms of the scaling properties of the k-th moments of the transition probability density 
p{x,t\O,0) for finding a configuration meta-particle centroid around position x at time t 
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given that it started at the origin at time zero. These are defined componentwise by 

m[^\t) = E [x^'{t)] ~ kel,2 (20) 
for the Xi coordinate, where the scahng exponent 6 determines the prevaihng type of trans- 



port at the time scale under consideration. The k-th moments (20) can be calculated from 
simulations of the individual-based model with the estimator 127 



w=^x:(x«(t))', (21) 

i=l 



where i = 1, . . . , Z is the number of simulation runs in the ensemble and each of the Xi{t) = 
Xi(0),Xi(l), . . . ,Xi(T) corresponds to a single time series of length T of group centroid 
positions along the xi coordinate observed at discrete time intervals of length At. 
The simplest possible scenario for collective-decision making occurs when the transfer of the 
bias of the informed sub-population leads on average to motion with effectively constant 
velocity Vi along the informed direction, which we identify without loss of generality with 



the positive xi axis. In this case the mean displacement {k = 1, in (20)) scales linearly with 



time after a transient determined by the characteristic time scale Tc, the time to consensus 

m^^\t) vit, (22) 

where the degree of consensus c can be defined by the ratio of the mean group velocity vi 
to the individual particle speed s 

c=^. (23) 

s 

Values of c close to one result from a distribution of individual particle orientations concen- 
trated around the informed orientation. On the other hand, if these tend to be distributed 
uniformly on (0,7r], one should expect comparatively smaller values of c, since the individ- 
ual particle velocities tend to cancel each other in this regime, which we associate to poor 
consensus. 

Figure [7] shows the behavior of the effective velocity vi versus the coupling constant u for a 
SPP swarm of ten particles and various proportions of informed individuals. We see that in 
all cases the degree of consensus increases as a power law of the coupling constant u with 
a slope that decreases as informed individuals are added. This means that if there is an 
optimum group velocity in some appropriately defined sense, it can be reached collectively 
by two different avenues. One is to have a small number of informed individuals at a high 
coupling constant, and the other is to have a large number of informed particles with low val- 
ues of the coupling constant. The power law dependence implies that the difference in values 
of the coupling constant for these two strategies can span orders of magnitude. Therefore, if 
the cost of recruiting informed individuals is less than that of leading the group, a possible 
optimal strategy, would consist of recruiting additional informed individuals, each of them 
with comparatively smaller values of u, instead of simply increasing the coupling constant 
of the informed group size, which comes at the additional complication of increasing the 
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time t (sees. ) 

Figure 6: Mean displacement m^^\t) along the informed direction xi versus time for a wide 
range of values of the coupling constant u. Results are averaged over 3 x 10^ independent 
simulation runs. Initial configurations are given by uniformly distributed locations within a 
circle of radius 0.5, and uniformly distributed orientations on the unit circle. In all cases, 
the mean displacement increases asymptotically linearly with time, indicating motion with 
a constant effective speed. The dip at early times for u = 0.007 (red dash-dot line, first from 
the bottom upwards) constitutes a signature of the transient; however a larger ensemble is 
required in this case is required to capture it accurately, since the values involved are much 
smaller than those present for higher values of u. However, for our purposes, the linear 
long time behavior is clearly shown. The simulation parameters are N = 10 individuals, 
A^^ = l,cr = 0.1 radians, Vat = 5.0 m, r^^ = 1.0 m, s = 1.0 m/sec. 
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probability of having the group spht apart. 

Naturally, these two strategies to reach the same target group velocity, are likely to have 
different accuracies. This can be more readily detected in measures of spread along the mean 



value, like the second moment. For this purpose we use the msd {k = 2, in (20)), which 
can also detect very efficiently the various types of collective behaviors, either transient or 
asymptotic, that contribute to macroscopic transport. For instance, the time to consen- 
sus can be detected sharply by a transition from linear (diffusion-dominated) or anomalous 
scaling to a quadratic one at the point in time where advection begins to dominate 

mf\t)r^vlt\ (24) 

Figure |8] shows the msd along the informed direction xi for a group of 10 individuals, one 
of them informed, and various values of the coupling constant. There is a transition from 
linear to quadratic scaling for non-negative values of w at a characteristic time scale Tc. As 
one increases the coupling constant, the informed sub-population becomes more efficient 
at transferring their bias to the whole group, signaled by an earlier time to consensus. 
The transient regime appears to be anomalous (supperdiffusive) in both naive {u = 0) and 
informed {u > 0) configurations. The anomalous transient can be better detected by looking 
at the scaling of the second order fluctuations, which requires removing the mean value in 



the definition of the moments (20) for k = 2 



(25) 



If a = 1 in (25) for some set of time scales, then the fluctuations behave as classical diffu- 
sion [221 ES] . Values of a different from one are dubbed anomalous and can be of two main 
types: sub-diffusive or 'trapped-diffusion' for < a < 1, and super-diffusive (sub-ballistic) 
or 'enhanced diffusion', if 1 < a < 2. These anomalous behaviors signal the presence of 
fluctuations that have persistent correlations in space, time, or both at macroscopically rel- 
evant scales [nni EH]- Slip/stick dynamics dominate the early time behavior of the second 



moment of the fluctuations along the Xi (Figure 10 ) and X2 (Figure |9]) coordinates, where 
the transport is clearly anomalous. There is a (sub-ballistic) super- diffusive transient that 
eventually decays to classical diffusion at a characteristic time scale r^. This is due to the 
alternation of bursts of advective flights (the slip phase) due to polarization of the orienta- 
tions, that are interrupted when the polarization is lost and the group moves much more 
slowly (stick). Since the msd eventually becomes diffusive, the temporal correlations in the 
mean velocity induced by the polarization eventually decay at a characteristic time scale 



Ta (shown in Figure 10), after which the fluctuations are diffusive, with identical diffusion 
coefficients along both coordinates. 

The role of the informed sub-population is more nuanced along the X2 coordinate (see Figure 
9|, in the sense that for higher values of the coupling constant, there is a clearly detectable 
sub-diffusive regime between the early time super-diffusive transient and the diffusive regime; 
this is due to reversals in velocity that are more marked along the X2 direction. Thus, the 
introduction of informed individuals at high coupling constants induces an anisotropy not 
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Figure 7: Dependence of tlie drift velocity vi along the informed direction xi on the coupling 
constant uj for various values of the number of informed individuals Np = 1 (black), 3 (blue) 
and 5 (red) for a swarm of iV = 10 individuals. In all cases, the drift increases as a power law 
of the coupling constant, but the exponent decreases as the number of informed individuals 
increases, due to the upper bound of the group velocity imposed by the individual particle 
speed. 
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only in the values of the diffusion coefficients but also in the manner in which the correlations 
decay. Whereas the the ADEM ([T]) allows for anisotropies in the diffusion coefficients along 
the two directions, the memory kernel M{t) must be identical along each direction. This pre- 
cludes the use of the ADEM ([T]) for predicting the behavior of the full 2-dimensional density, 
which would require a separate treatment along each coordinate. We can still however use it 
for the marginals along each direction in order to predict the mean squared displacements, 
and to extract the transport coefficients. 

Other macroscopic measures of interest are related to the precision of the decision-making 
process. This quantity is directly linked to the strength of the fluctuations along the informed 
direction, relative to those of a fully naive swarm of the same total population size. Higher 
precision in decision-making naturally implies a distribution of centroid positions along xi 
that is highly concentrated around the mean value. This is measured by the magnitude of 



the diffusivity along xi. Figure [TO] shows estimates of the second order fluctuations along the 
informed direction for various values of the number of informed individuals N13 = 0, 1, 3, 5, 7 
and a fixed coupling strength of u; = 0.5. All the fluctuations scale asymptotically linearly 
and thus they are dominated by classical diffusion, but with a diffusivity that decreases as 
the number of informed individuals iV^ increases. This agrees with the intuition that the 
precision of the decision-making process should improve as informed individuals are added 
to a group of fixed size. Figure |9] shows the mean squared displacement along the direction 
perpendicular to the informed orientation X2 for the same parameters in Figure |8j We see 
that the diffusivity along this coordinate is numerically indistinguishable between naive and 
informed configurations for the various values of the coupling constant used. This suggests 
a definition of decision-making precision given by the ratio of the diffusivities along the two 
coordinates, 

p = ^ < 1. (26) 

Values of p close to one indicate a spread of the meta-particle positions that is comparable 
along both coordinates, and thus poor precision in decision-making, and the opposite situ- 
ation occurs for values of p that are significantly less than one. 



3 Continous time random walks and the advection- 
difFusion equation with memory 

In this section we briefly review known results about continuous time random walk models 
(CTRW), that were originally introduced to describe the random motion of a particle on a 
disordered lattice (or medium). The main innovation of CTRW theory consisted in allowing 
the lattice spacing and updating times to become random variables themselves. It can be 
shown [m 130] that when the distribution of lattice spaces has finite moments of all orders, 
the evolution of the transition density for the location probability density of a CTRW is given 
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msd along Xi 




Figure 8: Mean squared displacement along the informed direction xi versus time for various 
values of the coupling constant u, for a total group size = 10 and one informed individual 
for non-zero values of u. There is an anomalous transport regime at early times evidenced by 
the scaling exponent a ~ 1.7 (25). This regime decays asymptotically to classical diffusion 
(linear scaling) in the naive configuration. When u becomes positive the transport becomes 
advective at a characteristic time scale Tc, signaled by the quadratic scaling of the msd, this 
time scale becomes shorter with increasing coupling strength. 
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Figure 9: Mean squared displacement along the X2 coordinate versus time. There is a 
superdiffusive transport regime at early times followed by a sub-diffusive phase that appears 
at relatively high values of the coupling constant (w > 0.3). The anomalous transient gives 
way to diffusive transport asymptotically in all cases with numerically indistinguishable 
diffusion coefficients. 
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Figure 10: Mean squared displacement with drift removal (25) along the informed direction 
Xi versus time for various values of the informed sub-population size and a fixed value 
{u = 0.5) of the couphng constant and total population size (A^ = 10). The super-diffusive 
transient eventually gives way to classical diffusion in all cases at a characteristic time scale 
Ta- Both the scaling exponent a and the asymptotic diffusivity decrease as the number 
of informed individuals increases, which signals an increase in the precision of the decision 
making for fixed values of the coupling constant, following the introduction of additional 
informed individuals while the total group size remains unchanged. 
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by the generalized master equation (GME) [301 EHl E21 EI] also called the advection-diffusion 
equation with memory (ADEM) when there is a drift [HI [3]. This generalization of the 
classical random walk leads to an evolution equation for the transition probability density 
of the particle position that is non-local in time, since the flux depends on a weighted time 
average over the full past. The weighting function is commonly called a 'memory function', 
and results from the wide range of transition rates originating from the spatial disorder. This 
approach has been successfully applied in models of anomalous diffusion [52], that typically 
require the memory to decay algebraically instead of exponentially. Approaches based on 
the CTRW and GME methods however are more general, and allow for any functional form 
of the memory, provided that it can be normalized. 

As will be discussed in Section |4| the memory function plays a fundamental role, since 
it encodes macroscopic transport coefficients of interest, together with their characteristic 
time scales. We will exploit these properties of memory functions in order to estimate the 
transport parameters associated with the various types of collective behaviors arising in SPP 
models of collective motion, with and without informed individuals. 

3.1 Continuous time random walks 

Continuous time random walks (CTRW) [HI |30l |3l| are a generalization of classical random 
walks [221 [56] where the jump size Ax and updating time At are allowed to become random 
variables. Sample paths are generated by drawing the jump size ^, and waiting time r from 
the joint probability density ip^C,, r). The elapsed time t„ for such a walker after n steps is, 

n 

i=i 

and the position x„(t), for a 2-dimensional walk, 

n 

The probability of observing a walker at position x at time t given that it started at the 
origin at time zero is, 

p(x, t) = 5(x) ^{t)+ [ r) p(x -^,t-T)d^ dT, (27) 

JM2 Jo 

where the survival function \['(t) is the cumulative of the waiting time marginal density of 
^(e,r) 

^{t) = l- [ [ ^{^,r)d^dT. (28) 

JR2 Jo 

For the particular situation were the jumps and waiting times are decoupled, the joint density 
iP{C,,t) can be rewritten as iP{t)X{^), where ipi^) is the distribution of waiting times, and 
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A(^) is the distribution of jumps. This assumption simphfies (27) to 

p(x,t)=(5(x)vl/(t)+ [ A(e) !\{T)p{-K-i,t~T)didT. 



(29) 



If the jump density A(^) has finite moments, and p(x, t) can be expanded in a Taylor series, 
it can be shown |3] that the differential version of (29) corresponds to an advection-diffusion 
equation generalized to non-local time. 



^P(x, t) 

dt 
p(x,0+) 



- f M {t - s) [va ■ Vp(x, s) - Da : V Vp(x, s)] ds 
Jo 

5(x), X G M^ t G M+ 



(30) 



where the memory term M{t) can alternatively be defined in terms of the Laplace transform 
of the waiting time density, or as the kernel of the velocity time autocorrelation (divided by 
2D so that it integrates to one). In the former case we have. 



M(6) 



(31) 



and /(e) denotes the Laplace transform of a function f{t) with e being the Laplace variable, 
and i is the characteristic time between transitions 



t 



Tip^r) dr. 



(32) 



The drift term va in (30) is related to the first moment of the jump pdf A(^), 

1 



Va 



t 



X A(x)(ix, 



(33) 



and the diffusivity tensor Da is given by the second moment of A(.^) 



D, 



2t 



xx"^ A(x)(ix. 



(34) 



In the drift-free case, the Laplace domain solution of the ADEM (30) is given by 

1 



p(x, e) 



27rM(e)/D7^ 



Kn 



M(e) 



X2 ,-y« 2 

1 _l_ -^2 



(35) 



which assumes that the off-diagonal components of the diffusivity tensor are zero, Di and 
D2 are the diffusivities along the xi and X2 coordinates, Kq is the modified Bessel function 
and M(e) is the Laplace transform of the memory with e being the Laplace variable. If the 
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drift vector has a single non-vanishing component which coincides with the Xi direction, the 
solution is [3] 



p(x, e) 



1 



27r M(e)/D^ 



X? + 



Do 



4- 



'M{e)vf 



(36) 



where Vi is the magnitude of the drift along the Xi coordinate. Finally, the Laplace domain 



expression for the mean squared displacement of the ADEM ( 30 ) along the direction of the 
drift is 

2vf —2/ N , 2£'i 



m 



Mtie) + —^Mi{e), 



(37) 



which yields all the characteristic time scales after Laplace inversion. 



4 Multiscale Method 

We start with the assumption that the meta-particle random walk follows an unknown 
CTRW with independent jump and waiting time distributions. In this case one can assert 
that after the velocity-autocorrelation equilibrates, the evolution of the transition density 
for the meta-particle location p{'x,t\O,0) is given by an advection-diffusion equation with 
memory (30) under relatively mild assumptions. The ADEM would be fully specified if an- 
alytical forms of the jump and waiting time densities were known on the basis of the SPP 
formulation. Unfortunately, this is not the case. We instead estimate them from a single 
velocity and centroid position time series obtained from a simulation run of the spp model 
with a combination of non-parametric and parametric methods. We show that this simple 
estimation procedure predicts mean squared displacements that are indistinguishable nu- 
merically from those estimated from an ensemble average over a large number of simulation 
runs. We exploit these results in order to explore a wide region of the parameter space, and 
obtain analytical results for the time to consensus based on the functional forms used in the 
parametric estimation of the memory. These results are exact for the case of exponential 



and Gamma density (40) memories, but only approximate for the truncated Mittag-LefHer 
case (|42|. 



4.1 Estimation of M{t) 

Although the memory in (pO) ) is defined in terms of the Laplace transform of the distri- 
bution of waiting times (31), a more convenient definition relates it to the time velocity 



autocorrelation of the random walker 



M{t) 



2{Di + D2, 



E[(v(r)-/x) 



V r 



t) - ^^)] 



(38) 
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where /i is the expected value of the random velocity v(t). The definition of the ADEM 



(30) allows different diffusivities along each coordinate, but not anisotropies in which the 



correlations decay differently along each component of the velocity. In general though, one 
should consider the memory separately along each component for an accurate description of 
the evolution of the transition pdf. Unfortunately, this turns out to be the case for informed 
swarms at high coupling constants, as can by seen after observing the differences for high 
values of the coupling constant in the mean squared displacements along X2 (Figure [9]) and 
the drift-corrected msd along xi (Figure 10). The former has a distinct sub-difussive regime 
that is not apparent in the latter. For the purposes of collective-decision making, it suffices 



to focus on the behavior of the msd (37) along the informed direction, and the macroscopic 
transport coefficients Di,D2 and Vi. 



We first compute a non-parametric estimate of the velocity auto-correlation function 
from a velocity time series {fi, f2, . . . , vt} obtained from a single simulation run of the SPP, 
where each of the Vi,i = 1, . . .T is the component of the meta-particle velocity along the 
informed direction, sampled at discrete time intervals At, and T is the length of the time 
series. We used the unbiased estimator |17] 



C{r) 



1 



T 



T 



i=l 



T-1, 



(39) 



where r is the time lag and v is the sample mean 

T 

T 



i=l 



The tabulated function that results from the non-parametric estimate (39) is fed to a non- 



linear least squares routine that yields a parametric estimate of the memory (see below). 
The Laplace transform of this function is then substituted into the expression for the mean 



squared displacement (37), or the transition pdfs (35) and (36), all of which can then be 



inverted numerically. The parametric estimate requires a 'template' function for the velocity 
auto-correlation that fits the data well and has a known analytical Laplace transform. We 
identified two functions that provide remarkably good fit and have very simple transforms. 



For lower values of the coupling constant this template is the Gamma density (Figure 11) 



fit) 



(40) 



r(i-/3) 

where r controls the exponential decay, and the exponent /3 controls the initial algebraic 



decay. The Laplace transform of (40) is simply 



1 

Ta 



(41) 
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The second function is appropriate for higher values of the couphng constant oj which leads 



to oscillations (see Figure 11). In this function the initial power law decay in the Gamma 



density in (40) is substituted by an exponentially truncated Mittag-LefHer function ^91159] . 

m = ^-^^ t^-' E^A-tVn) exp(-t/r.), (42) 

TeTa 

where the Mittag-LefHer function Ea^/3{—t'^/T^) is defined as 

-l)'=(r/r,)'^ 



E, 



a,l3\ 



k=0 



r{ak + f3) 



(43) 



where a and /3 are shape parameters and controls the transition between the early time 
and the asymptotic regime. The Laplace transform of the truncated Mittag-LefHer function 



(42) is also very simple 



-bt 



{b + e 



+ a 



(44) 



Figure 11 shows the results the estimation procedure for a swarm of = 10 individuals, one 



of them informed. The black marks show the non-parametric estimates based on (39) and 



the blue lines show the parametric fits using the Gamma density (40) for lower values of the 
coupling constant {uj = 0.1 and 0.3), and the truncated Mittag-LefHer function for higher 
values {u = 0.45 and 0.6). In all cases, the template functions provide remarkably good fit, 
including the oscillation that appears for higher values of u. The functions eventually decay 
to a constant value that corresponds to the drift squared vf, which we do not remove from 
the estimator, in order to be able to resolve the changes in the qualitative behavior of the 
correlation function for various values of u. Table [2] shows the parameter estimates in all 
four cases, together with goodness of fit values. Data for the velocity time series starts 
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SSE 


0.10 


0.0373 


0.012 




8.33 




0.20 


0.996 


2.3 X 10-6 


0.30 


0.0373 


0.038 




0.94 




0.24 


0.998 


2.9 X 10-^ 


0.45 


0.0326 


0.068 


0.55 


0.26 


1.09 


0.86 


0.999 


4.4 X 10-^ 


0.60 


0.0282 


0.094 


0.62 


0.37 


1.61 


0.79 


0.999 


2.3 X 10-6 



Table 2: Parameter estimates and goodness of fit values for the correlation functions in 
Figure 11 using the Gamma density (40) and truncated Mittag- LefHer function 



(42) 



as 



fitting templates. * The time scale is displayed in units of time for comparison with the 
exponential relaxation r^. SSE stands for Sum of the Squared Errors. 

being collected after a transient of 1000 time steps, after which time the time series becomes 



second-order stationary. Figure 12 shows that after a very short transient of a few hundred 
time steps, the estimators become very narrowly bounded and no trend with time is evident. 
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Figure 11: Estimated velocity autocorrelation function from a single time series of 1x10^ 
time steps (marks) and fitted functions (blue lines). The group size in the simulation was 
10 individuals, one of them informed. The black markers correspond to the non-parametric 
estimates, for various values of u. The continuous lines show the parametric fits with a 
Gamma kernel /(r) + vf 



g{T) + vf (|42| for u 
found in Table |2] 



40) for CO = 0.1 and 0.3, and a truncated Mittag-Leffler function 



0.45 and 0.6. Parameter estimates and goodness of fit values can be 
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Figure 12: Velocity autocorrelation function at zero lag Co(to) computed from a window of 
fixed length T = 1 x 10^ time steps, and shifting the origin of the first data point in the 
window To time steps from the absolute origin of the simulation run. Each graph corresponds 
to a different value of the coupling constant. The swarm simulation consisted of a total 
group size of = 10 individuals, of which one is informed. The arrow indicates the point at 
which data started to be collected for the estimates of the macroscopic transport parameters 
(to = 1000 time steps.) 
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4.2 Estimation of va and Da 

The drift coefficient vi can be estimated in two ways. The most straightforward is from the 
sample mean of the velocity time series {vi, V2, ■ ■ ■ ,vt}, 



vl = ^J2''^ (45) 



i=l 

and the other is based on the ffist moment of the jump kernel, 



f A = I / xi \i{xi)dxi (46) 



where i is the mean time between transitions (32), and Ai is the marginal of the jump 
kernel A(x) along the informed coordinate. Since the characteristic time t is not known, the 
estimator requires sampling the jump kernel A at various lags r. The characteristic time will 
be the value of r for which the estimator saturates, 

nr) 

where A(xi;r) is the sub-series of position differences along the direction xi sampled at 
time lag r from the position time series {xi,X2i ■ ■ ■ , xt}, where T is the total length of the 
series, and T(r) is the length of the sub-series sampled at lag r. Of course, the quality of 
the estimator decreases with r, because the length of each sub-series is twice as short as the 
preceding one. The lag dependent drift is then given by 

- f \ ^ 
^^a(^) = - 
r 

and the characteristic time can be calculated as the smallest value of the lag r* for which 
the equality 

that relates both estimators holds. Since a parametric form of the memory is already available 
(see Section 4.1), the diffusivity can be estimated from the Kubo-Green relationship [37 t l32 | 
[31] that relates transport parameters to time correlation functions, 

PCO 

D,= E [{v,{0) - V,) {v,{r) - V,)] dr. (48) 

^0 

It can also be estimated from the second moments of the jump kernel 

Di = ^ I {xi - fiif Ai(xi) dxi, (49) 
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where an estimator of Di is developed in a similar vein as that of the drift vx 



T{r) 



Di{r) = 2r(TiT)-l) ^ ' ' ^^^^ 



1=1 



the diffusion coefficient is the value for which -Di(t) reaches a plateau. The behavior of 



both estimators for the swarm meta-particle is shown in figure 14 The upper panel shows 
the results for naive configurations of various total population sizes, = 10 (red), = 50 
(green) and = 100 (blue). The dotted black line corresponds to the estimate of the 



diffusivity from the velocity time auto-correlation using the Kubo-Green relationship (48) 



and the rugged lines of various colors correspond the estimates of the diffusivity based on 



(50) that vary with the sampling lag r. The lower panel shows the comparisons between 
both methods for informed configurations of the same total population sizes as in the upper 
panel, but including informed individuals for the same coupling constants. In all the cases 
the proportion of informed individuals p = N/Np = 0.3 was kept constant. We observe that 
both methods converge to approximately the same value, in both naive and informed config- 
urations. We note that the characteristic time -the time at which the estimator saturates- 



increases with group size. The width of the oscillations in the estimator ( 50 ) increases with 
the lag r due to the finite size of the location time series, since for larger values of r, the 
number of data points used in the estimator decreases. 



4.3 Estimation of the time to consensus Tr 



We defined crudely the time to consensus as the time scale that determines the onset 
of the quadratic scaling in the mean squared displacement (Figure [s]) along the informed 
direction, which in the Laplace domain is given by 



m 



(51) 



where the coefficients vi and Di can be determined from (45) and (48) respectively. The 
parameters of the memory are calculated by the method described in Section 4J^ Given 
that the analytical Laplace transforms are known for both memory templates (40) and (42), 
substituting the Laplace transform of the Gamma memory (41) into (51) leads to 

2 



~ (2)/ \ 



(52) 



Likewise, for the substituting the Laplace transform of the truncated Mittag-Leffler function 
(|44|) yields 



my(e) 



2r, (r-i + e 



,a-2l3 



re (r-1 + e) 



(53) 
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Figure 13: Behavior of the estimator of the diffusion coefficient D{T) based on the Kubo- 
Green relationship (48) versus the length T of the meta-particle velocity time series. Open 
squares denote the value of the estimator using a truncated Mittag-Leffler kernel template 
for the velocity auto-correlation, and black circles correspond to a Gamma density. In both 
cases the dotted lines are the 95% confidence intervals. 
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Figure 14: Diffusion coefficients estimated via the Kubo-Green relationship (48) (dotted 
black lines), and from the variance of the jump kernel (50) sampled at various time lags 
r . Panel (a) shows estimates for purely naive swarms of total population size = 10 
(red), = 50 (green) and iV = 100 (blue). Panel (b) shows the estimates for informed 
configurations. The three cases share the same fraction of informed individuals p = 0.3 and 
coupling constant uj = 0.3, the total population sizes are color coded as in panel (a). We 
note that both methods succeed in providing the asymptotic value of D. There is an overall 
reduction in diffusivity as the total group size increases. The diffusivity also decreases in 
informed groups compared with naive ones of the same total size. 
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Figure 15: Comparison between the mean squared displacement along the informed direction 
Xi (20) estimated from an ensemble of 3000 simulation runs (black marks) and that obtained 
from the inverse Laplace transform of the msd (37) based on the fitted ADEM (blue con- 
tinuous lines), with parameters estimated from a single simulation run. We used a Gamma 
density memory kernel for the lower values of the coupling constant {u = 0,uj = 0.1) and an 
exponentially truncated Mittag-LefHer function for the remainder cases {u = 0.3, u = 0.6). 
In all cases the total population size consisted of = 10 individuals, and informed configu- 
rations consisted of one informed individual in all cases. 
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In the case of the msd for the Gamma density memory (52 ), it is possible to invert analytically 



the Laplace transform. The msd with the truncated Mittag-Leffler memory can be inverted 
numerically using the inversion algorithm of de Hoog [18j. Before we can use the results of 



the analytical and numerical inversions of (52) and (|53|) we show in Figure 15 comparisons 



between the msd obtained from an ensemble of simulation runs of the swarm meta-particle 
(black marks) and that obtained by inversion of the Laplace transforms of the mean squared 



displacements (52) and (53) (blue lines) based on the ADEM assumption, with parameters 



estimated from a single simulation run of the SPP, using the method outlined in Sections 4.2 



and 4A In all cases the method based in the ADEM is able to capture accurately both the 
transient and the asymptotic behavior. In order to use these results to calculate the time 
to consensus, Tc we first note that in the simpler case of a memory of the form of a Dirac 



distribution 6{t), Laplace inversion of (51) is straightforward 



m 



in which case Tc is the smallest time scale for which the contribution due to advection is 
larger than that of diffusion, 

vl > 2Di t, 

which leads to 

r. = ^. (54) 



A similar procedure can be carried out for non-trivial choices for the memory. The first of 
these is an exponential memory with a relaxation time scale Tq = 1/6. This functional form 
dominates the asymptotic behavior in both the Gamma and the truncated Mittag-Leffler 
memory kernels if the anomalous time scale in the latter is sufficiently fast compared with 
1/b. An analogous procedure yields the time to consensus 



2Di (l 



l-Cxp(-bTe) 
bTc 



bTc 



+ 



exp(-6rc 



(55) 



which requires an iterative solution. The full Gamma kernel (40) results in 
2D, (l + ^ + ^ [exp(-6r,) - (/3 + br, - 1) S^{bT,)] 



6+4br,(/3-l)+2/3(2/3-5) , (bTe)2(i-/^) 



+ 



r(4-2/3) 



[exp(-6re)(2/3 + br, - 1) - R{b, /3, rj] 



where 



R{b, P, rj = (6 + + 46r,(/3 - 1)) + 2/3(2/3 - 5) ^2/3-3(K) 
and Sa{x) is the exponential integral 



(56) 



(57) 



£Jx) 



oo g-xt 



1 



dt. 



38 



Unfortunately, we were unable to find an analytical inversion of the Laplace transform of 



(53) for a Mittag-LefHer memory kernel. However, both memory kernels are dominated 



asymptotically by the exponential truncation. For simplicity, we used the exponential ap- 



proximation (55) of the time to consensus for the macroscopic analysis of the efficiency 
of collective decision making for various group sizes, values of the coupling constants and 
proportions of informed individuals. 

4.4 Results 



Figure [T6| shows estimates of the three key macroscopic parameters of swarm meta-particles. 
The magnitude of the diffusivity Di along the informed direction (left column), the drift 
Vi (center column), and the time to consensus (right column, logarithmic scale) for three 
total population sizes = 10 (top row), iV = 50 (center row) and = 100 (lower row). In 
all the graphs the horizontal axis corresponds to the coupling constant u G [0, 0.6], and the 
vertical axis to the relative fraction p of the informed population size to the whole group. 
We see that the precision of the collective decision, measured by the ratio of the diffusivities 



along both coordinates ( 26 ) increases with the coupling constant and the number of informed 



individuals. Similarly, the degree of consensus (23), measured by the ratio of the drift vi to 
the individual particle speed, increases as well with the coupling constant and the informed 
fraction. Smaller groups move faster than larger ones, but at the cost of a loss in precision. 
Finite size effects are of paramount importance in this class of problems. Given that the 
diffusivity decreases with group size as was also detected before traditional approaches 
where macroscopic quantities are calculated in the limit of very large population sizes are 
not particularly useful in this context. The time to consensus Tc decreases with increasing 
number of informed individuals and coupling strength. This is not surprising since it is tied 
to first order to the ratio Di/vf. Interestingly, it appears to be invariant to group size and 
controlled by the time scale of the exponential relaxation which increases as the group 
size grows. 



5 Final comments 

This study suggests that both the transient and the asymptotic regimes of swarming popu- 
lations -with strong alignment and in the presence of an orientation bias- can be concisely 
approximated by an advection-diffusion equation with memory. The presence of an orien- 
tation bias together with macroscopic bursts of alignment, alternating with an unpolarized 
phase, lead to quite non-trivial time correlations in the mean group velocity, which persist 
over macroscopically relevant time scales. These must be explicitly accounted for in order to 
capture accurately the macroscopic parameters that typify the various collective states to- 
gether with their characteristic time scales. This observation is consistent with recent results 
by Griinbaum et al [25] who found that local-in-time advection-diffusion equations even with 
density dependent coefficients could not fully capture the fluxes of individual-based models 
of swarming populations when alignment was an important contributor to the dynamics at 
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Figure 16: Estimates of the diffusion coefficient (left column) along the informed direction D, 
the mean group speed (center column) f , and time to consensus (right) r for various values 
of the proportion of informed individuals p (vertical axis), coupling constant cu (horizontal 
axis), and total population sizes. The first row {D\,v\,r\) corresponds to the case N = 10, 
the second {D2, V2, T2) to = 50 and the third {D^, v^, t^) to N = 100. 
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the level of the individual particle. That study focused on looking at the fluxes of fission- 
fusion populations, without informed individuals, for various values of the density in order 
to try to find a functional form that fitted the dependence of the transport coefficients on 
the population density. We explored a much more limited range of population sizes, but 
instead looked in more detail at the temporal dependence of the mean squared displacement, 
and the various transport behaviors shown at each time scale. Of course, both methods are 
not in opposition but complement each other. In the future, we would like to integrate both 
approaches in such a way that both the density-dependence and memory effects are included 
in a single transport model of swarming populations with alignment. 

We find that the mean group velocity increases as a power law of the couphng constant, and 
that the exponent of the power law decreases as the number of informed individuals increases. 
We also find -in agreement with earlier work [T3]- that the total group size has a dramatic 
impact in the collective transport properties. Smaller groups tend to move with higher ve- 
locities, but at the expense of a higher diffusivity and thus less precise decisions. This may 
have important implications for evolutionary studies of simple models of collective-decision 
making, where there is presumably costs associated with recruiting informed individuals into 
the population, by having a relatively high value of the coupling constant and by making 
erroneous decisions (Vishwesha Guttal et al, personal communication). If some value of 
the mean group velocity along the informed direction is optimal in a way that maximizes a 
measure of individual-level fitness, there are a number of possible ways to achieve it. One 
possible path is to have a small number of informed individuals, each with a relatively high 
coupling strength, while another is to have a larger number of informed individuals but with 
a much smaller coupling strength. A very rich trade-off space is likely to occur in this class 
of systems, particularly if one allows for variability in total population size. 
Remarkably, the efficiency of collective decision-making, understood as the time scale at 
which an effective drift becomes detectable over the diffusive component of the meta-particle 
random walk, seems to be invariant with respect to group size. What seems to determine 
the efficiency is a combination of the fraction of informed individuals and the strength of 
the orientation bias. This arises from the fact that this quantity ultimately depends on the 
ratio D/v"^ and the characteristic time scale of the exponential decay in the memory (55). 
The time velocity auto-correlation emerges from the ADEM approach as the key macro- 
scopic summary statistic. It quantifies the relative contributions to macroscopic transport 
from each collective behavior, and allows the specification of their characteristic time scales. 
Although the ability of time correlation functions to connect microscopic dynamics with 
observed macroscopic regimes has been known in non-equilibrium statistical physics for at 
least four decades since the seminal work of Kubo [36] , Mori |13] , Green [21] , Zwanzig [61] , 
MontroU [H] and Kenkre ^D], to our knowledge it is a relatively unexplored concept in 
movement and spatial ecology, where Markovian models have dominated the scene 0B] , per- 
haps with the notable exception of correlated random walks [TOl ESj Il5] . We would like to 
emphasize a subtle point though, which is that the temporal memory of the ADEM does 
not necessarily imply that the individual walker has information about the past in order 
to make movement decisions about the future. The memory arises naturally as a result of 
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the ensemble average of a continuous time random walk in the presence of a wide range of 
transition rates. These can result from internal properties -like an updating clock with a 'fat 
tail' instead of an exponential one- or external factors such as behavioral variability due to 
complicated social interactions or spatial structure in the landscape that results in slip/stick 
dynamics; these can occur quite naturally if there are corridors with preferential directions 
of motion alternating with regions where movement can be described with Brownian motion. 
We believe that this ecological interpretation of the time velocity auto-correlation function 
is likely to be useful not only to unravel the connections between individual-based models of 
movement and dispersal and their continuum approximations as we have seen in this study, 
but also for other areas of ecology where interdependencies between an individual organ- 
ism's dispersal strategy, spatial heterogeneity in the landscape, and temporal variability in 
resource availability become intertwined in observed individual trajectories, particularly in 
the nascent field of movement ecology. 

Future work will be devoted to a generalization of the SPP model to density-dependent 
asynchronous updating, in the sense that each of the social interactions is associated with 
an exponential clock that is parameterized by the local density, in a similar way to what is 
done in locally regulated models of plant population dynamics with spatial structure O E] . 
Given that the CTRW-ADEM can predict the the full density and not just the first moments, 
future work will be devoted to this issue in order to explore first passage times. We will also 
explore the situation when there are two conflicting preferential directions, where it remains 
to be seen whether the ADEM has the capability of capturing the bifurcations that have 
been detected in individual-based simulations [T3]. This will require generalizations of the 
ADEM involving anisotropy in the memory kernel. 

6 Acknowledgements 

The authors are grateful for the support received from the National Science Foundation 
(Award ID:EF-0434319) and DARPA (Award ID:HR001-05- 1-0057). Insightful discussions 
with I. D. Couzin, L.Giuggioli, V.M. Kenkre and F.Bartumeus are gratefully acknowledged. 



References 

[1] I. Aoki. A simulation study on the schooling mechanism in fish. Bulletin of the Japanese 
Society of Scientific Fisheries, 48(8):1081-1088, 1982. 

[2] Igor S. Aranson, Audrey Sokolov, John O. Kessler, and Raymond E. Goldstein. Model 
for dynamical coherence in thin films of self propelled microorganisms. Physical Review 
E, 75:040901(R), 2007. 



42 



[3] Brian Berkowitz, Andrea Cortis, Marco Dentz, and Harvey Scher. Modeling non-Fickian 
transport in geological formations as a continuous time random walk. Reviews of Geo- 
physics, 44(2): Art. No. RG2003/2006, 2006. 

[4] Eric Bertin, Michel Droz, and Guillaume Gregoire. Boltzmann and hydrodynamic de- 
scription for self-propelled particles. Physical Review E (Statistical, Nonlinear, and Soft 
Matter Physics), 74(2):022101, 2006. 

[5] Benjamin Bolker and Stephen W. Pacala. Using moment equations to understand 
stochastically driven spatial pattern formation in ecological systems. Theoretical Popu- 
lation Biology, 52:179-197, 1997. 

[6] Benjamin Bolker, Pacala Stephen W., and Simon A. Levin. Moment methods for eco- 
logical processes in continuous space. In Ulf Dieckmann, Richard Law, and Johan A.J 
Metz, editors. The Geometry of Ecological Interactions, volume 1 of Gamhridge Studies 
in Adaptive Dynamics, pages 388-411. Cambridge University Press, 2000. 

[7] J Buhl, DJT Sumptcr, ID Couzin, JJ Hale, E Dcspland, ER Miller, and SJ Simpson. 
From disorder to order in marching locusts. Science, 312(5778):1402-1406, JUN 2 2006. 

[8] Andrea Cavagna, Alessio Cimarelli, Irene Giardina, Alberto Orlandi, Giorgio Parisi, 
Andrea Procaccini, RafFaele Santagati, and Fabio Stefanini. New statistical tools for 
analyzing the structure of animal groups. Mathematical Biosciences, 214(1-2) :32-37, 
JUL-AUG 2008. 

[9] Yao-Li Chuang, Maria R. D'Orsogna, Daniel Marthalcr, Andrea L. Bertozzi, and Lin- 
coln S. Chayes. State transitions and the continuum limit for a 2D interacting, self- 
propelled particle system. Physica D, 232:33-47, 2007. 

[10] Edward A. Codling, Michael J. Plank, and Simon Benhamou. Random walk models in 
biology. Journal of the Royal Society Interface, 5(25):813-834, AUG 6 2008. 

[11] R. R. Coifman, l.G. Kcvrckidis, S. Lafon, M. Maggioni, and B. Nadlcr. Diffusion 
maps, reduction coordinates, and low dimensional representation of stochastic systems. 
Multiscale Modeling & Simulation, 7(2):842-864, 2008. 

[12] Larissa Conradt and Christian List. Group decisions in humans and animals: a survey 
introduction. Philosophical Transactions of the Royal Society B-Biological Sciences, 
364(1518):719-742, MAR 27 2009. 

[13] Larissa Conradt and Timothy J. Roper. Group decision-making in animals. Nature, 
421(8):155-158, 2003. 

[14] Andrea Cortis, Claudio Gallo, Harvey Scher, and Brian Berkowitz. Numerical simulation 
of non-fickian transport in geological formations with multiple scale heterogeneities. 
Water Resources Research, 40:W04209, 2004. 



43 



[15] Ian D. Couzin, Jens Krause, Nigel R. Pranks, and Simon A. Levin. Effective leadership 
and decision-making in animal groups on the move. Nature, 433:513-516, 2005. 

[16] Ian D. Couzin, Jens Krause, Richard James, Graeme D. Ruxton, and Nigel R. Pranks. 
Collective memory and spatial sorting in animal groups. Journal of Theoretical Biology, 
218:1-11, 2002. 

[17] Andras Czirok, Albert-Laszlo Barabasi, and Tamas Vicsek. Collective motion of self- 
propelled particles: Kinetic phase transition in one dimension. Physical Review Letters, 
82(1):209-212, 1999. 

[18] P.R. de Hoog, J.H. Knight, and A.N Stokes. An improved method for numerical in- 
version of Laplace transforms. SIAM Journal on Scientific and Statistical Computing, 
3(3):357-366, 1982. 

[19] R. Eftimie, G. de Vries, and M. A. Lewis. Complex spatial group patterns result from 
different animal communication mechanisms. Proceedings of the National Academy of 
Sciences of the United States of America, 104(17) :6974-6979, APR 24 2007. 

[20] Radek Erban, loannis G. Kevrekidis, and Hans G. Othmer. An equation-free computa- 
tional approach for extracting population-level behavior from individual-based models 
of biological dispersal. Physica D- Nonlinear phenomena, 215(l):l-24, MAR 1 2006. 

[21] G Plierl, D Griinbaum, S Levin, and D Olson. Prom individuals to aggregations: The 
interplay between behavior and physics. Journal of Theoretical Biology, 196(4) :397-454, 
PEB 21 1999. 

[22] C.W. Gardiner. Handbook of Stochastic Methods. Springer- Verlag, second edition, 1985. 

[23] Wayne M. Getz and David Saltz. A framework for generating and analyzing movement 
paths on ecological landscapes. Proceedings of the National Academy of Sciences of the 
United States of America, 105(49):19066-19071, DEC 9 2008. 

[24] MS Green. Comment on a paper of mori on time-correlation expressions for transport 
properties. Physical Review, 119(3) :829-830, 1960. 

[25] Daniel Griinbaum, Karen Chan, Elizabeth Tobin, and Michael T. Nishizaki. Non- 
linear advection-diffusion equations approximate swarming but not schooling popula- 
tions. Mathematical Biosciences, 214(1-2) :38-48, JUL- AUG 2008. 

[26] Daniel Griinbaum and Akira Okubo. Modeling social animal aggregations. In Simon A. 
Levin, editor. Frontiers in Theoretical Biology, volume 100 of Lecture Notes in Biomath- 
ematics, pages 296-325. Springer- Verlag, 1994. 

[27] Paul R. Halmos. The theory of unbiased estimation. Annals of Mathematical Statistics, 
17(l):34-43, 1946. 



44 



[28] JW Haus and KW Kehr. Generalized effective-medium approximation for particle - 
transport. Physical Review B, 36(10): 5639-5642, OCT 1 1987. 

[29] A. Huth and C. Wissel. The simulation of the movement of fish schools. Journal of 
Theoretical Biology, 156:365-385, 1992. 

[30] V. M. Kenkre. The generalized master equation and its applications. In Uzi Landman, 
editor. Statistical Mechanics and statistical methods in theory and application, pages 
441-461. Plenum Pubhshing Corporation, 1977. 

[31] V. M. Kenkre, Z. Kalay, and P. E. Parris. Extensions of effective-medium theory of 
transport in disordered systems. Physical Review E, 79(1, Part 1), JAN 2009. 

[32] V.M. Kenkre. Memory formalism, nonlinear techniques and kinetic equation approaches. 
In V.M. Kenkre and K. Lindenberg, editors. Modern Challenges in Statistical Mechan- 
ics: Patterns, Noise and the Interplay of Nonlinearity and Complexity, pages 63-102. 
American Institute of Physics, 2003. 

[33] loannis G. Kevrekidis, C. WiUiam Gear, James M. Hyman, Panagiotis G. Kevrekidis, 
Olof Runborg, and Constantinos Theodoropulos. Equation-free, coarse-grained mul- 
tiscale computation: enabling microscopic simulators to perform system-level analysis. 
Communications in Mathematical Sciences, 1(4): 715-762, 2003. 

[34] J Klafter and R Silbey. Derivation of the continuous time random walk equation. 
Physical Review Letters, 44(2):55-58, 1980. 

[35] AUison Kolpas, Jeff Moehlis, and loannis G. Kevrekidis. Coarse-grained analysis of 
stochasticity-induced switching between collective motion states. Proceedings of the 
National Academy of Sciences, 104(14):5931-5935, 2007. 

[36] R Kubo. Statistical-mechanical theory of irreversible processes .1. general theory and 
simple applications to magnetic and conduction problems. Journal of the Physical 
Society of Japan, 12(6):570-586, 1957. 

[37] Ryogo Kubo, Morikazu Toda, , and Natsuki Hashitsume. Statistical physics II: Nonequi- 
librium Statistical Mechanics. Springer- Verlag, 1991. 

[38] A. Merrifield, Mary R. Myerscough, and N. Weber. Statistical tests for analysing di- 
rected movement of self-organising animal groups. Mathematical Biosciences, 203(1) :64- 
78, SEP 2006. 

[39] R Mctzlcr and J Klafter. The restaurant at the end of the random walk: recent devel- 
opments in the description of anomalous transport by fractional dynamics. Journal of 
Physics A-Mathematical and General, 37(31) :R161-R208, AUG 6 2004. 

[40] A Mogilner and L Edelstein-Keshet. A non-local model for a swarm. Journal of Math- 
ematical Biology, 38(6):534-570, JUN 1999. 



45 



[41] E.W. MontroU and H. Scher. Random walks on lattices iv. continuous time random 
walks and influence of absorbing boundaries. Journal of Statistical Physics, 9(2), 1973. 

[42] Sung Joon Moon, B. Nabet, Naomi E. Leonard, Simon A. Levin, and I. G. Kevrekidis. 
Heterogeneous animal group models and their group-level alignment dynamics: An 
equation-free approach. Journal of Theoretical Biology, 246(1):100-112, MAY 7 2007. 

[43] H Mori. Statistical-mechanical theory of transport in fluids. Physical Review, 
112(6):1829-1842, 1958. 

[44] Benjamin Nabet, Naomi E. Leonard, Iain D. Couzin, and Simon A. Levin. Dynamics of 
Decision Making in Animal Group Motion. Journal of Nonlinear Science, 19(4):399-435, 
AUG 2009. 

[45] Ran Nathan, Wayne M. Getz, Eloy Revilla, Marcel Holyoak, Ronen Kadmon, David 
Saltz, and Peter E. Smouse. A movement ecology paradigm for unifying organismal 
movement research. Proceedings of the National Academy of Sciences of the United 
States of America, 105(49):19052-19059, DEC 9 2008. 

[46] Akira Okubo and Simon Levin. Difussion and Ecological Problems: Modern Perspec- 
tives, volume 14 of Interdisciplinary Applied Mathematics. Springer- Verlag, second edi- 
tion, 2001. 

[47] S.J. Orfanidis. Optimum Signal Processing. An Introduction. Prentice-Hall, second 
edition, 1996. 

[48] Juha K. Parrish, Steven V. Viscido, and Daniel Gunbaum. Self-organized fish schools: 
An examination of emergent properties. Biological Bulletin, 202:296-305, 2002. 

[49] Igor Podlubny. Fractional Differential Equations. Mathematics in Science and Engi- 
neering. Academic Press, 1999. 

[50] Craig W. Reynolds. Flocks, herds and schools: A distributed behavioral model. Com- 
puter Graphics, 21 (4): 25-34, 1987. 

[51] Audrey Sokolov, Igor S. Aranson, John O. Kessler, and Raymond E. Goldstein. Concen- 
tration dependence of the collective dynamics of swimming bacteria. Physical Review 
Letters, 98:15802, 2007. 

[52] David J. T. Sumpter. The principles of collective animal behavior. Philosophical trans- 
actions of the royal society B : Biological Sciences, 361:5-22, 2005. 

[53] David J. T. Sumpter, Jens Krause, Richard James, Iain D. Couzin, and Ashley J. W. 
Ward. Consensus decision making by fish. Current Biology, 18(22): 1773-1777, NOV 25 
2008. 



46 



[54] John Toner, Yuhai Tu, and Sriram Ramaswamy. Hydrodynamics and phases of flocks. 
Annals of Physics, 318(1) :170-244, 2005. 

[55] Chad M. Topaz, Andrea L. Bertozzi, and Mark A. Lewis. A nonlocal continuum model 
for biological aggregation. Bulletin of Mathematical Biology, 68(7):1601-1623, OCT 
2006. 

[56] N.G. Van Kampen. Stochastic Processes in Physics and Chemistry. North-Holland, 
2001. 

[57] Ashley J. W. Ward, David J. T. Sumpter, Lain D. Couzin, Paul J. B. Hart, and Jens 
Krause. Quorum decision-making facilitates information transfer in fish shoals. Proceed- 
ings of the National Academy of Sciences of the United States of America, 105(19):6948- 
6953, MAY 13 2008. 

[58] B.J. West, P. Grigolini, R. Metzler, and T.F. Nonnenmacher. Fractional diffusion and 
levy stable processes. Physical Review E, 55(1, Part A):99-106, JAN 1997. 

[59] Bruce J. West, Mauro Bologna, and Paolo Grigolini. Physics of fractal operators. 
Springer, 2003. 

[60] Christian A. Yates, Radck Erban, Carlos Escudero, Iain D. Couzin, Jerome Buhl, loan- 
nis G. Kevrekidis, Philip K. Maini, and David J.T. Sumpter. United by noise: Ran- 
domness helps swarms stay together. Proceedings of the National Academy of Sciences, 
106(14) :5464-5469, 2009. 

[61] Robert Zwanzig. Time-correlation functions and transport coefficients in statistical 
mechanics. Annual Review of Physical Chemistry, 16:67-&, 1965. 



47 



